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Abstract. A model for homogeneous anisotropic incompressible turbulence is proposed. The model generalizes the GISS 
model of homogeneous isotropic turbulence; the generalization involves the solution of the GISS equations along a set of 
integration paths in wavenumber (k-) space. In order to make the problem tractable, these integration paths ("cascade 
^s^j ' lines" ) must be chosen in such a way that the behaviour of the energy spectral function along different cascade lines should 

be reasonably similar. In practice this is realized by defining the cascade lines as the streamlines of a cascade flow; in the 
simplest case the source of this flow may be identified with the source function of the turbulence. Owing to the different 
approximations involved, the resulting energy spectral function is not exact but is expected to give good approximative 
values for the bulk quantities characterising the turbulent medium, and for the measure of the anisotropy itself in particular. 

The model is then applied to the case of low Prandtl number thermal convection. The energy spectral function and the 
bulk quantities characterizing the flow are derived for different values of the parameter S = Ra<r. The most important new 
finding is that unlike the anisotropy of the most unstable mode in linear stability analysis the anisotropy of the turbulence 
does not grow indefinitely with increasing S but it rather saturates to a relatively moderate finite asymptotic value. 
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1. INTRODUCTION 



In recent years it has become clear that a better understanding of the anisotropy of turbulence in stellar 
convective zones is important for several reasons. It is now widely believed that the anisotropy of 
turbulence plays a key role in governing the transport of angular momentum in the convective zones of 
stars and therefore in determining their differential rotation (Riidiger, 1989). More accurate knowledge 
of the anisotropy could therefore constrain the theories of differential rotation, which at present involve a 
large number of free parameters. Furthermore, it is possible to show (Petrovay, 1990) that the intriguing 
morphological properties of solar and stellar convection found in recent numerical simulations (Stein 
and Nordlund, 1989, Nordlund and Dravins, 1990) can be predicted on the basis of sufficiently general 
convection theories if the anisotropy is known. From its influence on differential rotation and on the 
morphology of convection it follows that the anisotropy of turbulence is important in understanding the 
dynamo mechanism (Brandenburg et al. , 1991, Petrovay, 1991). Finally, some current theories of stellar 
convection (e. g. Xiong, 1980) rely on the assumption that the anisotropy of turbulence remains moderate 
throughout the convective layer, while other formalisms (Canuto 1989, 1990) include the anisotropy as a 
free parameter to be specified. All this suggests that a method to compute the anisotropy of turbulence 
that is simple enough to be applicable in practice while it still retains the essential physical ingredients 
of the problem could be very useful in several astrophysical areas. 

It seems to be plausible to treat homogeneous anisotropic turbulence as a first step. Despite the fact 
that homogeneous anisotropic turbulence is "to some extent a non-problem" (Leslie 1973), in the sense 
that any process generating (anisotropic) turbulence relics on the presence of an inhomogeneity on the 
large scales while the turbulence becomes isotropic on the small scales, it has been widely recognized that 
the study of this subject is of high importance for several reasons. First, it is expected that with increasing 
k wavenumber the turbulence tends to homogeneity faster than to isotropy (Leslie 1973), so there exists 
a wavenumber-regime where the model applies. Secondly, in certain flows (e. g. in grid turbulence 
sufficiently far from the grid) the homogeneous anisotropic approximation may be valid even on the 
scale of the energy-containing eddies. Thirdly, and perhaps most importantly, the study of homogeneous 
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Fig. 1. The anisotropy of the "typical" (T) and the strongest (S) mode may be similar or totally different, depending on 
how sharply peaked the spectral function is. 



and anisotropic turbulence may be an important intermediate step towards the fully general case of 
inhomogcneous anisotropic turbulence. 

As efforts to tackle the problem of anisotropic turbulence on the basis of two-point closure models like 
the DIA did not bear much fruit owing to the intimidating complexity of the equations, Canuto et al. 
(1990) propose that, to begin with, simpler heuristic models should be used, followed by the construction 
of more complex models (in analogy with the historical development of the theory of homogeneous 
isotropic turbulence). 

The measure of anisotropy is usually characterised by quantities like 

or 

s = =§ (2) 

v 

where l x , l y , l z are "characteristic scales" of the turbulence (e. g. macroscales or correlation lengths) 
in different directions, v is the turbulent velocity and z is a preferred direction (e. g. due to gravity). 
In contemporary practice values of these quantities are often "guessed" on the basis of their values for 
the most unstable mode in linear stability analysis. Turbulent values for the diffusivities arc sometimes 
used in the analysis as an attempt to take the nonlinear interaction of modes into account. It is however 
rather uncertain whether the anisotropy computed in this way represents the anisotropy of the turbulence 
correctly, not only because nonlinear mode interactions may be expected to distort the shape of the 
spectrum shifting the maximum, but also because quantities like (1) or (2) are defined as averages over 
all wavenumbers, weighted e. g. by 

<Mk) = J exp(-zke)^(r)« J (r + 0d 3 e, (3) 

the velocity covariance spectral function (or its trace F = $m), and not as values for the single mode 
corresponding to the maximum of F (k). The difference is highlighted in Figure 1 (for the sake of clarity 
in one dimension). Depending on how sharply peaked F(k) is, the "typical" k value may be similar to 
or totally different from the location of the maximum of F(k). 

In fact, for the case of low Prandtl number convective turbulence in numerical experiments the value 
of the s parameter defined in (2) is found to be very nearly constant at a moderate value of s ~ 0.37 for 
a range of the 

S = (4) 
X 2 



ANISOTROPY OF TURBULENT CONVECTION 



3 



parameter covering several orders of magnitude (Chan and Sofia, 1989). (Ra is the Raylcigh number, a 
the Prandtl number, g the gravity acceleration, S the volume expansion coefficient, (3 the superadiabatic 
temperature gradient, x the thermodiffusion coefficient^], D the layer depth, usually identified with the 
Hp pressure scale height.) This is in contrast with the anisotropy of the most unstable mode in linear 
stability analysis which grows indefinitely with increasing S (Canuto, 1990). On the other hand, Canuto 
(1990) has shown that if molecular diffusivities are replaced by turbulent ones in the linear analysis, 
the anisotropy of the most unstable mode, though still quite high and slowly increasing with 5, will be 
seriously reduced. This was the first theoretical result indicating that nonlinear interactions of modes 
have the potential to limit the anisotropy of turbulence. 

The present paper proposes a method to derive the measure of the anisotropy as well as other bulk 
quantities in a consistent (though approximate) way. Section 2 contains the general description of the 
method. This model is then applied to the astrophysically relevant case of low Prandtl number convective 
turbulence in Section 3. The results, summarized in Section 4, agree with the above mentioned numerical 
simulations. Section 5 concludes the paper. 

2. THE MODEL 

2.1. General description of the model 

The model presented below is a generalization of the GISS model developed by Canuto, Goldman and 
Chasnov (1987, hereafter CGC). The GISS model is a single-point closure scheme belonging to the class 
of "turbulent viscosity approach" models; it has been shown to compare well with experiments and DIA 
results (see CGC). 

We suppose that <&jj(k) is a "smooth" function, i. e. that the spectrum of the allowed modes is not 
discrete. In cases of astrophysical interest, for which this model was primarily developed, this is often 
a reasonable approximation as the studied turbulent regions here form parts of much larger systems (e. 
g. stars) and so the eigenmodes populate the k-space densely enough. The assumption may also be 
reasonably valid in some laboratory experiments (e. g. grid turbulence). However, for systems with a 
discrete and sparsely populated mode spectrum the present approach cannot be applied. 

In its principal axis frame the velocity covariance spectral function takes the form 

/$i \ 

$ = $ 2 . (5) 
V $ 3 / 

For each $i(k) component we have an rts^k) source function or growth rate, supposed to be known from 
the linear stability analysis. (The dissipation rate is also supposed to be built into ns,i-) So we can form 
an hs source tensor function in such a way that in the principal axis frame of $ the form of ns is 

(n s ,i \ 
n S = n s ,2 . (6) 
V nsfl) 

The total energy input per unit time into a given k mode is then n<j(k) : $(k). (A : B = AijBij.) 

The basic equation of the model expresses the equality of the energy fed into a given V volume of the 
k-space with the energy transferred away from that volume: 

f n s (k') : ®(k')d 3 k' = f d 3 k' [_ T(k',k")d 3 k" , (7) 
JV Jv Jv 

V being the complement of V. The T(k',k ) transfer term is to be determined by the closure. In the 
turbulent viscosity approach the r. h. s. of equation (6) is taken to be the product of the total vorticity 
of modes in V with the turbulent viscosity. 

* \ = \F m \(pc p \\7T\)~ 1 with p the density, c p the specific heat at constant pressure, T the temperature and F m the 
microscopic (radiative + conductive) energy flux. 
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In the case of isotropy equation (7) for a continuous manifold of spheres centered on the origin of the 
coordinate frame uniquely determines 3>, once a = 4>/F is specified. In the case of anisotropy, however, 
we need several (in general three) such intersecting manifolds of surfaces. Equivalcntly, one may take 
a continuous manifold of curves originating from a common P starting point and covering the k-space. 
This manifold defines a "curvilinear polar coordinate frame" with coordinates k, 9c, <\>c ( k is the length 
along each curve measured from P, the 9c and 4>c "angles" effectively index the curves). Now in equation 
(7) we may choose for V a "curvilinear cone" with angle Afic centered around the < n' < k stretch of 
the (9c, 4>c) curve. Taking the ASlc — ► limit and using the notation k 2 as defined by 



we get 



where 



Jo 



d 3 k = k 2 dKd(l c (8) 
h s {n\9c,(j>c) ■ §(k' ,9c,4>c)k 2 c (n' ,9c,(j>c)dn' = y(K,9c,<f>c)vt(K,9c,<i>c) (9) 
y(n 1 9 c ,4>c)= I F( K ',9 c ,^c)k 2 (K\dc,<f>c)k 2 c ^',0c,<f>c)d K ' (10) 



is the total vorticity along the (0, k) stretch of the curve. 

According to the turbulent viscosity concept v t can be written in the form of an integral of independent 
contributions. On a purely dimensional basis 

Vt {K,9 c Ac)= I ^J-d 3 k", (11) (11) 

Jv(k,0 c ,<Pc) n c (\s. ) 

n^ 1 being a characteristic timescale of the turbulence. The GISS model closure is 

in c = k 2 v t (12) 

(7 = (2/3Ko) 3 where Ko is the Kolmogorov constant). 

Equation (11) shows that in order to solve the system of equations along a given path, the value of F 
must be known along all the other paths. How can we avoid this awkward recursion? 

Suppose that our choice of the integration paths is so fortunate that the behaviour of F(k) = 
F(k 7 9c,4>c) is very similar (or even identical) for all (9 c Ac) paths. In the case of this "fortunate" 
choice we will call the paths cascade lines. Equation (11) can now be written as 

v t {K,9cAc) = ** r k2 c( K "> c, 4>c) F( /V t°\ dn", (13) 
Jk n c (n ,9 C , (pc) 

i. e. the integration over the angles could be replaced by a 47r factor and the recursion has disappeared. 
The method of finding the cascade lines (i. e. such a fortunate choice of paths) will be discussed in the 
next subsection. 

Now the equations (9), (10), (12) and (13) can be transformed into a more convenient form, following 
the treatment in CGC with slight modifications. 

Differentiating (9) with respect to k and using the notation a — ^/F we get 

(hs ■ a) + Airy-ric 1 = v t k 2 . (14) 

From this 

2 7 n c = {h s : a) + [(h s : a) 2 + 16^ 7 y] 1/2 . (15) 
Still following CGC we introduce 

V = Any + §7«c-(16) (16) 

With this, equation (15) becomes 

7 n c = \(h s : a) + [±(n s : a) 2 + | 7 V] 1/2 . (17) 
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Now differentiating (16) and using the closure (12) 

dV 2 1 dk 2 

^=^T 2 ^ (18) 

(Partial derivatives are taken with (8c, 4>c) =const. , i. e. along the cascade lines.) 

Were a(K,8c,4>c) an d k 2 (n,6c,4>c) known, equations (16)-(18) would determine y(n, 8c,4>c), an d 
therefore F(k, 9c, <j>c)- On the basis of the original integral form of the equations the initial condition is 

V(k = 0) = ^-(n s :&f. (19) 

For a and k 2 , the asymptotic behaviour in the limits k — ► and k — > oo is known (for incompressible 
turbulence). In the intermediate regimes interpolation formulas may be used. These, and the choice of 
the cascade lines will be discussed in the next subsection. 



2.2. Choice of the cascade lines and quenching 

In the isotropic case the cascade lines may obviously be chosen as straight lines starting from the origin 
(Figure 2a). The behaviour of F(k) is in this case exactly identical for the different cascade lines. In 
the anisotropic case, however, this choice is not appropriate anymore. If, for instance, Tr(ng : a) has a 
sharp peak in a certain regime of the k-space, marked by a dashed outline in Figure 2 (a quite common 
case in physical processes generating turbulence) the energy input along paths crossing this regime will 
greatly exceed the input along other paths, so the F(n) values will be very different for different paths at 
large k. On the other hand, one might expect that a choice of the cascade lines that follows the general 
structure of the ng(k) function defining the problem (like the one sketched in Figure 2b) will lead to a 
much more balanced distribution of the energy input between different paths. 

These considerations lead us to specify the cascade lines the following way. Let the lines be the 
streamlines of a j c vector field. As for large k-s we expect isotropy, and therefore radial cascade lines, we 
can put rot j c = or j c = Vip c . The divergence of j c is now equated to the source function: 

V 2 Vc = -Tr(n s ). (20) 

An equally defensible approach would be to put 

\7 2 yj c = -(h s : a)F (21) 

on the basis that it is the total energy input that determines the shape of F(k). As however no formula 
can in general ensure that F(n) will be exactly identical for every path and equation (22) would involve 
much more numerical work (as using it coupled to the equations (16)— (18) some kind of iterative solution 
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would be necessary for F), it is more economic to use equation (20). Nevertheless, alternative choices of 
the cascade lines can be used to check the dependence of the results on this particular choice. 

Physically we expect that with increasing k wavenumber the turbulence becomes more and more 
isotropic. This property enables us to find the asymptotic form of k 2 (k) and d(k) for large k values. 

For k 2 we simply expect k 2 — > k 2 as k — > oo; while at points in the immediate neighbourhood of P 
we expect k 2 — > k 2 = (k — kp) 2 for truly 3-dimensional anisotropy and k 2 — > k 2 Q = nk Xt p\k — kp| for 
axially symmetric turbulence (i. e. with only one preferred direction). In the intermediate regime we may 
simply use an arbitrary interpolation formula describing the "quenching" of the anisotropy towards large 
k-s: 

k 2 c = -^k 2 + -^—k 2 c0 (22) 
l+p l+p c '° ' 

where 

Q is a numerical factor which must be order of unity as kp is basically the single important "length scale" 
in k-space. 

In order to be able to find a similar formula for a we must restrict ourselves to incompressible flows 
for which $ has only two nonzero independent components, say $1 = (cf. Batchelor, 1953). In this case 
a will take the form 



&= I 1-a I (24) 
0a 

and in the k — > 00 limit a — ► 1/2 while in the neighbourhood of P where (in the framework of the present 
model) all the energy in the modes comes directly from the energy input by the instability generating 
the turbulence we have a = ns,s/{ns,2 + n s,3)- Our interpolation formula will be 

a = 0.5- Kt ■ r- . (25) 

l+p (ns,2 + ns,3) 1 +P 

The precise functional form of these interpolation formulas is physically not too relevant, the sole 
important parameter that influences the values of the bulk properties is the Q quenching parameter. 

2.3. Discussion of the model 

To summarize what was said above, the recipe to compute F(k) for homogeneous anisotropic incompress- 
ible turbulence goes as follows. 

Starting from the fis source function defining the problem one solves the 

V 2 Vc - -Tr(ns) (26) 

Poisson-equation for the ip c cascade potential; physically we expect ip c — > e/k as k — > 00, so the boundary 
condition at large k values may be chosen as 

^ e(M , + ^^=0. (27) 

Now taking the gradient of ip c one can compute a desired number of cascade lines starting from the P 
maximum site of tp c . The system 

dV 2 1 dk 2 

^=^W (28) 

7 n c = l(h s : a) + [\{h s : a) 2 + §7^] 1/2 . (29) 
V = Any + \ in 2 c . (30) 

F = ihjk (31) 



ANISOTROPY OF TURBULENT CONVECTION 



7 




Fig. 3. Owing to the assumption of a preferred direction of energy transport, a spurious narrow "hole" appears near the 
maximum of -F(k) in the present model. 



can now be integrated along each cascade line with the initial condition 

V( K = 0) = ^(n s :a) 2 . 



(32) 



The expressions for fc 2 and a are given by equations (22)-(25). 

The resulting F(k) can then be used to compute the integrals defining the bulk quantities charac- 
terising the turbulence, among them the anisotropy parameters such as (1) or (2). Let us derive e. g. 
the expression for the s parameter defined in equation (2). It is a straightforward algebraic exercise to 
compute the components of a in the fixed k x k y k z frame of reference (in the incompressible case, using 
the kQcj) polar coordinates): 



a xx — a (cos 2 9 cos 2 (j) — sin 2 (f>) + sin 2 
a xy — a [(cos 2 9 + 1) sin cos 0] — sin cos < 
ol X z — cz sin 9 cos 9 cos (f> 
a yy = a (cos 2 9 sin 2 — cos 2 <f)) + cos 2 <f) 



(33) 
(34) 



a sin 



Using this, for s we get 



2s + 1 



Jo°° dr Jo # Jo <*(r, 0, <f>)F{r, 6, 0) sin 9 d9 
Jo°° dr $T W Jo F ( r > e > 4) sin3 9 d0 



(35) 



It is obvious that the resulting F(k) is only a rather crude approximation of the real velocity covariance 
spectral function. An assumption inherent in applying equation (6) to describe the flow is that the energy 
cascade has a preferred direction; this will lead to a spurious "hole" in F(k) around P where F(k) goes 
to zero (Figure 3). This "hole" is narrow enough for its effect on the bulk quantities not to be serious; 
it basically corresponds to the lack of backscatter in the isotropic GISS model. On the other hand, the 
fact that some cascade lines may actually run towards decreasing k values (cf. Figure 2b) means that the 
anisotropic model does contain some backscatter — an improvement with respect to the isotropic case! 

The approximation in equation (12) has the effect that the asymptotically isotropic behaviour of F(k) 
at large k values is not strictly guaranteed anymore. As the bulk properties of the turbulence are mainly 
determined by the energy-containing eddies at low k values, this inconsistency is not expected to cause 
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serious errors in the bulk quantities. (Yet a quenching of the resulting F(k) function with the same Q 
parameter as in equation (22) may be used to correct this flaw.) 

Altogether the method described here seems to be appropriate for a first approximative derivation 
of the values of the bulk quantities, most notably the anisotropy itself, in a homogeneous, anisotropic, 
incompressible turbulent flow. 



3. APPLICATION 

For the numerical implementation of the method we nondimcnsionalize every quantity by choosing D = 
X= 1: 

k = kD n s — n s D 2 x~ 1 

y = yD 2 x ~ 2 k= kD (36) 
F = FD^x^ 2 V = VD 4 X ~ 2 

etc. In order to simplify the notation, the tildes will be omitted from here on, but we must remember 
that all the formulae below refer to the nondimensional quantities. 

The problem at hand is axially symmetric, so in x, z, <j) cylindrical coordinates we may restrict ourselves 
to the k x k z -p\a,ne and n s has the form n s , xx = n s , zz = n s ; ns,^ = n s , x< p = n s , z <j, = n s . xz = 0. The 
formula for ns(k) was written down in the CGC paper; for low Prandtl numbers it is 



n s (k) = -h 2 



1 h 2 

4 k 2 



(37) 



As in practice we would like to apply this formula to stratified layers (identifying D = Hp), a k co cutoff 
wavenumber is introduced so that for values \k z \ < k co , ns is taken to be zero. This is intended to imitate 
the decrease of the growth rate with increasing vertical eddy size above a size of ~ 2Hp in linear stability 
analysis (Narasimha and Antia 1982; see also the Figure 1 in Chan and Serizawa 1991). In order to check 
the dependence of the results on the choice of k co , the computations were performed with two different 
k co values (1.0 and 3.14). 

The a(k) and fc^(k) functions are given in equations (22) and (25); the Q parameter values studied 
were 0.5, 1, 2, 3 and oo. 

The equations of the model were solved on a 32 x 32 rectangular grid in a square-shaped regime of 
the k x k z plane: < k x , k z < k max = 100. As the most interesting behavior of F(k) is expected at low k 
values, the grid was not uniformly spaced; instead, new variables 

k x = \og(k x + 1) 
k y = \og{k y + 1) 

k = log(« + 1) (38) 
k = log(fc + 1) 

were introduced and the grid spacing was chosen to be uniform in these variables. In the new variables 
equations (27)-(31) and (26) take the form 

dV 2 1 dk 2 

m=^Y 2 ^ (39) 

7 n c = \(n s a) + [i(n s «) 2 + § 7 V] 1/2 (40) 
V = ATTy+\ in 2 c (41) 

F-i'-Lg? (42) 
k z k£ exp k ok 

V(k = 0) = ^(n s a) 2 (43) 

cxp(-2^)^ + cxp( _ 2 fc z )^ - CX p(-2^)^ - cxp(-2fc z )|^ = ~2n s , (44) 
dk 1 . Okj. Ok x Ok z 
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Fig. 4. Contour lines of the -F(k) double energy spectral function on the k x k z plane for the case S = 10 10 , k co = 1, Q = 1. 



The boundary conditions for the ip c cascade potential are 

cit kx — 
at k z = 

at k x = 100 ^ 
at k z = 100. 

The elliptic problem defined by (43) and (44) was solved for ij) c by a 7-point multigrid iteration scheme 
using Numcrcal Algorithms (NAG) software library routines. After computing the resulting cascade lines, 
the (38)-(42) GISS equations were integrated along each cascade line with the Euler-Cauchy method. 
The resulting F values were then used to interpolate F on the mesh points of a regular rectangular and 
a polar grid and the bulk quantities characterizing the turbulent flow (weighted integrals of F(k)) were 
computed. 

The whole procedure was performed for many different values of S in the regime 10 2 -10 17 (the range 
of S in the solar convective zone) and with the k co and Q values quoted above. In order to check the 
dependence of the results on the choice of the cascade lines, in some cases alternative definitions for ip c 
were also used (as described in Section 2.2). As expected, the results proved to be fairly insensitive to 
the choice of the cascade potential as long as it by and large reflected the general geometrical structure 
of the source function. 



ft=° 
^ = 

(fcmax+l)fc max , S^c _ n 

[(cxpfc x -l) 2 + (cxpfcz-l) 2 ]i/2 Vc "T" d'k z 

(fcmax+l)fcrriax / , dfg _ Q 

[(0X P fc x -l)2 + (0 X pfc z -l)2]l/2 ^ C Qk x 



4. RESULTS 

The contour lines of F(k) are shown for an example run in Figure 4. The preference of high k x /k z modes 
over low k x /k z modes (i. e. vertical anisotropy) is apparent, particularly at low k values (in the range of 
the energy-containing eddies). Figure 5 shows the cascade lines for the same example. As explained in 
Section 2.4., the present model does not ensure the isotropic behavior of the solution as k — > oo, but this 
has little effect on the bulk properties. 
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Fig. 5. Cascade lines and contour lines of the cascade potential (dashed) for the same case as in Figure 4. 



logF 




Fig. 6. The F(k) angle-integrated double energy spectral function for the case 5 = 10 , k co = 3.14, Q = 1. Note the 
backscattcr at low wavenumbers. 



In Figure 6 the 



F(k) = 2n [ F(k,6)k 2 sin6 d6 (46) 
Jo 

integrated double energy spectral function is shown. This figure is to be compared with the Figure 23 
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log S 



Fig. 7. The s anisotropy parameter as a function of S for different Q values. 



of the CGC paper. Contrary to the isotropic GISS model, some backscatter is noticeable in the present 
model, though it is not so emphasized as in the DIA. 

Values of the s anisotropy parameter for different models are shown in Figure 7 as a function of S. 
These values can be computed from 

_ dk f* ' a{k,e,<t>)F{k,e,4>)smedid 

2s ~~\~ 1 — — 



f* max dk F(k, 6, (f>) sin 3 6 d6 



(47) 



i. e. the 2D version of Equation (34). It is apparent that instead of going to zero, s saturates to a finite 
limiting value as S — * oo. This qualitative behavior is independent of the choice of the Q quenching 
parameter, though the actual asymptotic value of s does depend on Q. The value s = 0.37 found in 
Chan and Sofia's (1989) simulations corresponds to Q ~ 3-4 but the anisotropy remains quite moderate 
even with Q — > oo. A similar effect, with a weaker dependence on Q is observed in the shape of the x(S) 
function (Figure 8). The x quantity defined in (1) was calculated here using 

C i = 2 7 r/ dk I F{k)k T k z sm.6 d6 



rkmax rir 

1~ 1 = 2tt dk F(k)k x k 2 sm6( 

Jo Jo 

Z" 1 = 2tt / dk F(k)k z k 2 sine d9. 
Jo Jo 



10 Jo 

In Figure 9 the 5-dependence of the nondimensional convective flux 



dk I n s (k)a(k)F(k)k z sinede 
'o 



(48) 



(49) 



is shown. For comparison, the CGC result and some mixingdength theory (mlt) calibrations (with the 
choice I = Hp = D as length scale) are also shown. It should be noted though that a direct comparison 
of the present Boussinesq model with mlt models of the strongly stratified solar convective zone is not 
necessarily sensible. The computed curve lies near the lower limit derived by Spruit (1974) and is by a 
factor of 7 lower (at large S values) than the isotropic CGC curve. Similarly, the total specific kinetic 
energy 



K = 7T f 

Jo 

(Figure 10) and the turbulent viscosity 

v t = 2tt / 
Jo 



dk f F(k)k 2 sinede 
Jo 

F(k)k 2 



dk 



I 

Jo 



n c (k) 



sin 6d6 



(50) 



(51) 
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Fig. 8. The x anisotropy parameter as a function of S in the present model (continuous line) and for the most unstable 
linear mode on the basis of the formulas of Canuto (1990) (dashed line). 




- 



2 4 6 8 10 12 14 16 

logs 

Fig. 9. The $ nondimcnsional convective flux as a function of S for the present work, for the CGC model, for the mlt 
formula of Gough and Weiss (1976), and the observational limits deduced by Spruit (1974) (dashed area). The mlt fluxes 
were nondimcnsionalizcd assuming I = Hp. 



(Figure 11) arc also considerably lower than the isotropic results. From these results the coefficient of 
the K-e relation 

e = &K 2 /v t (52) 

is found to be £3 = 0.184, which is a factor of two higher than the experimental value. Similarly, for 
Smagorinsky's constant we get C ~ 0.001, to be compared with the latest experimental value of 0.003 
(quoted by CGC). As in CGC, the agreement of the computed coefficients with the experimental ones to 
within a factor of 2 or 3 indicates that the qualitative behavior of ns(k) is typical for a general class of 
mechanisms generating turbulence. 
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Fig. 10. The K specific turbulent kinetic energy as a function of S. 
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Fig. 11. The ut turbulent viscosity as a function of S. 

5. CONCLUSION 

The application of the homogeneous anisotropic incompressible turbulence model desribed in Section 2 
has shown that the anisotropy of turbulent thermal convection does not grow indefinitely with incrcsing 
S, but it rather saturates to a finite (and moderate) limiting value. This result is in agreement with 
numerical experiments but it is in contradiction with earlier theoretical expectations based on linear 
stability analysis. It remains to be seen how the inhomogeneity (nonlocality) of turbulence will influence 
these results, but it appears that the assumption of near isotropy is not the main factor limiting the 
validity of present-day convection theories. 

The computations show that the resulting <!> convective flux is practically insensitive to the value of 
the Q quenching parameter, unless Q <C 1. The physical reason for this is that the flux is predomi- 
nantly determined by the low wavenumber energy-containing eddies, so details of the quenching at higher 
wavenumbers do not have serious effects on it. The comparison with numerical experiments shows that Q 
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is indeed of order unity. The low sensitivity of $ to Q means that the ^'optoS 1 ^ 2 scaling remains valid for 
any form of the (unknown) Q(S) function. As the Q parameter basically determines the efficiency of the 
coupling between different velocity components (the "toroidal" and "poloidal" components of Massaguer, 
1991), it follows that the secondary instabilities generating toroidal motions do not significantly influence 
the energy transport. 

The work reported here has been the first attempt to study anisotropic turbulent convection in a 
self-consistent (though approximate) way. It is needless to say that further, more sophisticated models 
are necessary to extend these studies, in particular to the inhomogeneous case. 
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